function[v,y] = generate_data(n,mb,mc,G,whether_randomize,whether_write_data)

% function for generating data for the beta estimating approach
% n = number of individuals
% mb = number of binary explanatory variables
% mc = number of continuous explanatory variables
% G = function used to generate data (typically inverse link function, 
%      but it's possible to use something different from the actual model 
%      which is solved using Newton-Raphson)

% v = n X (mb+mc) set of explanatory variable data
% y = binary outcome data

% Set the random seed (so we get reproducible results)
reset(RandStream.getGlobalStream,1);
rand_init = 125677;
%rand_init = 863;
randn(1,rand_init);

m = mb+mc;
% Create a beta vector
beta_generator = randn(1,m); % use randn to get from a N(0,1) distribution
% Generate synthetic data
for i = 1: n
    vbinary = randi([0,1],mb,1);
    vcontinuous = rand(mc,1); % use rand -- uniform from [0,1]  
                              % The idea is that these are usually positive
                              % variables, and so can be scaled by a
                              % maximum value to convert to [0,1].
    vtemp = cat(1,vbinary,vcontinuous);
    ytemp(i) = G(beta_generator * vtemp);
    v(i,:) = vtemp';
end
if strcmp(whether_randomize,'y') % can make data highly random
                                 %---and hence not amenable to any regression 
                                 % model---by choosing a large noise amount
    noise_amount = 0.01;
    ytemp = ytemp + noise_amount * randn(1,n);
end
y = (ytemp>mean(ytemp))';
% randomise further by flipping a fraction
if strcmp(whether_randomize,'y')  % an additional randomization possibility
    flip_fraction = 0.02;
    flip_number = round(n * flip_fraction);
    flip_indices = randi([1 n],flip_number,1);
    flip_indicator = zeros(1,n);
    flip_indicator(flip_indices) = 1;
    flip_indicator = flip_indicator';
    y = flip_indicator .* (1-y) + (1-flip_indicator) .* y;
end


% Write data to file (if requested to do so)
if strcmp(whether_write_data,'y')
    % prepare formatting to save data to file
    filename = strcat("b-mb",num2str(mb),"mc",num2str(mc),"n",num2str(n),".txt");
    fileID = fopen(filename,'w');
    fmt_binary = repmat(' %2.0f',1,mb);
    fmt_continuous = repmat(' %12.8f',1,mc);
    fmt_binary_string = repmat(' %2s',1,mb);
    fmt_continuous_string = repmat(' %12s',1,mc);
    for mm = 1: m
        v_string(mm) = strcat("v",num2str(mm));
    end
    % write variable names line to file
    fprintf(fileID,' %2s','y');
    if mc > 1
        fprintf(fileID,fmt_binary_string,v_string(1:mb));
        fprintf(fileID,fmt_continuous_string,v_string(mb+1:end-1));
        fprintf(fileID,' %12s\n',v_string(end));
    elseif mc==1
        fprintf(fileID,fmt_binary_string,v_string(1:mb));
        fprintf(fileID,' %12s\n',v_string(end));
    else
        fprintf(fileID,fmt_binary_string,v_string(1:mb-1));
        fprintf(fileID,' %2s\n',v_string(mb));
    end
    % write data to file, a line at a time
    for i = 1: n
        fprintf(fileID,' %2.0f',y(i));
        if mc > 1
            fprintf(fileID,fmt_binary,v(i,1:mb));
            fprintf(fileID,fmt_continuous,v(i,mb+1:end-1));
            fprintf(fileID,' %12.8f\n',v(i,end));
        elseif mc == 1
            fprintf(fileID,fmt_binary,v(i,1:mb));
            fprintf(fileID,' %12.8f\n',v(i,end));
        else
            fprintf(fileID,fmt_binary,v(i,1:mb-1));
            fprintf(fileID,' %2.0f\n',v(i,mb));
        end
    end
    fclose(fileID);
end
